New method for the extrapolation of finite-size data to infinite volume 
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We present a simple and powerful method for extrapolating finite-volume Monte Carlo data to infinite volume, 
based on finite-size-scaling theory. We discuss carefully its systematic and statistical errors, and we illustrate it 
using three examples: the two-dimensional three-state Potts antiferromagnet on the square lattice, and the two- 
dimensional 0(3) and 0(oo) fr-models. In favorable cases it is possible to obtain reliable extrapolations (errors 
of a few percent) even when the correlation length is 1000 times larger than the lattice. 



Quantum field theorists are interested primar- 
ily in infinite systems; but Monte Carlo simula- 
tions must perforce be carried out on lattices of 
finite linear size L, limited by computer memory 
and speed. This raises the problem of extrapo- 
lating finite- volume data to _L = cxd. We present 
here a simple and powerful method for perform- 
ing this extrapolation, based on finite-size-scaling 
theory [1]; and we discuss carefully its systematic 
and statistical errors. We illustrate the method 
using three examples: the two-dimensional three- 
state Potts antiferromagnet on the square lattice 
[2], and the two-dimensional 0(3) and O(cxd) a- 
models [3,4]. We have found — much to our sur- 
prise — that in favorable cases it is possible to 
obtain reliable extrapolations (errors of a few per- 
cent) at ^/ L as large as 10-1000. More details can 
be found in [5]. 

Consider, for starters, a model controlled by 
a renormalization-group (RG) fixed point having 
one relevant operator. Let us work on a peri- 
odic lattice of linear size L. Let ^(p,L) be a 
suitably defined finite-volume correlation length 
(we use the second-moment correlation length de- 
fined by equations (4.11)-(4.13) of [6]), and let O 
be any long-distance observable (e.g. the corre- 



lation length or the susceptibility). Then finite- 
size-scaling theory [1] predicts that 

where fo is a universal function and w is a 
correction-to-scaling exponent. Hence, if s is any 
fixed scale factor (usually we take s = 2), 

ffi^ = Fa{aP,L)/L) + 0{r'',L--) (2) 

where Fq can be expressed in terms of fo, f^- 

Our method proceeds as follows [7]: Make 
Monte Carlo runs at numerous pairs (/3, L) and 
(j3,sL). Plot 0(j3,sL)/0(j3,L) versus C(j3,L)/L, 
using those points satisfying both ^(/3, L) > some 
value ^min and L > some value Lmin ■ If all 
these points fall with good accuracy on a sin- 
gle curve — thus verifying the Ansatz (2) for 
■f > Cmin , L > Lmin — choose a smooth fitting 
function Fo- Then, using the functions F^ and 
Fo, extrapolate the pair (C,0) successively from 
L ^ sL ^ s'^ L -^ . . . ^ cxD. 

We have chosen to use functions Fo of the form 



Fo{x) = 1 + aie-i/^ + ... + a^e' 
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(3) 
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This form is partially motivated by theory, which 
tells us that F{x) -^ 1 exponentially fast as a; ^ 



[10]. Typically a lit ol order 3 < n < 12 is 
sufficient; we increase n until the x^ of the fit 
becomes essentially constant. The resulting x^ 
value provides a check on the systematic errors 
arising from corrections to scaling and/or from 
the inadequacies of the form (3). 

The statistical error on the extrapolated value 
of Ooo(P) = 0(P,(X)) comes from three sources: 
(i) error on 0{p,L), which gets multiplicatively 
propagated to Ooo ; (ii) error on ^(p,L), which 
affects the argument x = ^(/3, L)/L of the scaling 
functions F^ and Fq ', and (iii) statistical error in 
our estimate of the coefficients ai, . . . ,a„ in F^ 
and Fo ■ The errors of type (i) and (ii) depend on 
the statistics available at the single point (p,L), 
while the error of type (iii) depends on the statis- 
tics in the whole set of runs. Errors (i)-|-(ii) [resp. 
(i)-|-(ii)-|-(iii)] can be quantified by performing a 
Monte Carlo experiment in which the input data 
at (P,L) [resp. the whole set of input data] are 
varied randomly within their error bars and then 
extrapolated. 

The discrepancies between the extrapolated 
values from different lattice sizes at the same /3 
— to the extent that these exceed the estimated 
statistical errors — indicate the presence of sys- 
tematic errors and thus the necessity of increasing 
Lmin and/or ^rnin aud/or n. 

A figure of (de)merit of the method is the rel- 
ative variance on the extrapolated value Ooo(P), 
multiplied by the computer time needed to obtain 
it. We expect this relative variance-time product 
[for errors (i)-|-(ii) only] to scale as 



RYTF (p,L) ^Cco(P) 



d+z. 



GoUoo(P)/L (4) 



where d is the spatial dimension and Zi„t,o is the 
dynamic critical exponent of the Monte Carlo al- 
gorithm being used; here Go is a combination 
of several static and dynamic finite-size-scaling 
functions, and depends both on the observable O 
and on the algorithm but not on the scale factor 
s. As ^oo/L tends to zero, we expect Go to di- 
verge as {(,oolL)~'^ (it is wasteful to use a lattice 
L ^ S,co)- As ^oo/L tends to infinity, we expect 
Go ~ (Coo /Ly [5], but the power p can be ei- 
ther positive or negative. If p > 0, there is an 
optimum value of^oo/L; this determines the best 



lattice size at which to perform runs for a given 
l3. If p < 0, it is most efficient to use the smallest 
lattice size for which the corrections to scaling are 
negligible compared to the statistical errors. 

Our first example [2] is the two-dimensional 
three-state Potts antiferromagnet on the square 
lattice, which is believed to have a critical point 
at l3 = cxD [11]. We used the Wang-Swendsen- 
Kotecky cluster algorithm [12], which appears to 
have no critical slowing-down (r.v,* a^2 < 5 

uniformly in /3 and L) [2]. We ran on lattices 
L = 32, 64, 128, 256, 512, 1024, 1536 at 153 difl^er- 
ent pairs {j3,L) in the range 5 < ^oo ^ 20000. 
Each run was between 2 x 10® and 2.2 x 10^ 
iterations, and the total CPU time was modest 
by our standards (about 2 years on an IBM RS- 
6000/370). We took S^rnin = 10 and Lmin = 128 
and used a quintic fit in (3); the result for F^ 
is shown in [2,5] (x^ = 75.41, 66 DF, level = 
20%). The extrapolated values from different lat- 
tice sizes at the same /3 agree within the esti- 
mated statistical errors (x^ = 43.03, 75 DF, level 
> 99%). The result for G^ is shown in [5]: the er- 
rors are roughly constant for ^oo/L > 0.4 but rise 
sharply for smaller ^oo/L. In practice we were 
able to obtain ^oo to an accuracy of about 1% 
(resp. 2%, 3%, 5%) at ^oo ~ 1000 (resp. 2000, 
5000, 10000). 

Next let us consider [3,4] the two-dimensional 
0(3) (T-model (see Caracciolo's talk for more 
details). We used the Wolff embedding algo- 
rithm with standard Swendsen-Wang updates; 
again critical slowing-down appears to be com- 
pletely eliminated. We ran on lattices L = 
32,48,64,96,128,192,256,384,512 at 180 difl^er- 
ent pairs {j3,L) in the range 20 < ^oo ~ 10®. 
Each run was between 10® and 5 x 10® iterations, 
and the total CPU time was 7 years on an IBM 
RS-6000/370. We took S^rnin = 20 and used a 
tenth-order fit. There appear to be weak correc- 
tions to scaling (of order < 1.5%) in the region 
0.3 < Cl/L < 0.7 for lattices with L < 64-96. We 
therefore chose Lmin = 128 for ^l/L < 0.7, and 
Lmin = 64 for Cl/L > 0.7. The result for F^ is 
shown in [4,5] (x^ = 72.91, 73 DF, level = 48%). 
The result for G^ is shown in [5]; at large ^oo/L it 
decreases sharply, with a power p Ri —2 in agree- 



ment with theory [5]. In practice we obtained ^oo 
to an accuracy of about 0.2% (resp. 0.7%, 1.1%, 
1.6%) at ^oo ~ 10^ (resp. 10^, 10^ 10^). 

We also carried out a "simulated Monte Carlo" 
experiment for the 0(N) cr-model at N = co, 
by generating data from the exact finite-volume 
solution plus random noise of 0.1% for L = 
64,96,128, 0.2% for L = 192,256 and 0.5% for 
L = 384,512 [which is the order of magnitude 
we attain in practice for 0(3)]. We considered 
35 values of /3 in the range 20 < ^oo ^ 10'^. 
We used ^min = 20 and Lmin = 64 (in fact 
much smaller values could have been used, as 
corrections to scaling are here very small) and 
a ninth-order fit; for two different data sets we 
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114 (resp. 118) with 166 DF. In prac- 



tice we obtain ^oo with an accuracy of 0.6% (resp. 
1.2%, 2%, 3%) at ^oo ~ 10^ (resp. 10^ 10^, 10«). 
Here we can also compare the extrapolated values 
C^^''(J3) with the exact values ^S?'"=*(/3). Defining 
^ = j:pr^"'(P)-C^'"HP)?/'T'(P), we find for 
the two data sets TZ = 17.19 (resp. 25.81) with 
35 DF. Only 6 (resp. 9) points differ from the ex- 
act value more than one standard deviation, and 
none by more than two. 

Details on all of these models will be reported 
separately [2,4]. 

The method is easily generalized to a model 
controlled by an RG fixed point having k relevant 
operators. It suffices to choose k — I dimension- 
less ratios of long-distance observables, call them 
R = (Ri, . . . , Rk-i); then the function Fq will 
depend parametrically on R(p,L). In practice 
one can divide i?-space into "slices" within which 
Fo is empirically constant within error bars, and 
perform the fit (3) within each slice. We have 
used this approach to study the mixed isovec- 
tor/isotensor cr-model, taking R to be the ratio 
of isovector to isotensor correlation length [3,4]. 

The method can also be applied to extrapolate 
the exponential correlation length (inverse mass 
gap). For this purpose one must work in a system 
of size L'^-^ xT with T > ^expiP, L) (cf. [8]). 
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